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Macroscopic  models  of  radiative  transfer  as 
applied  to  computation  of  the  radiation  field  in 
the  solar  atmosphere 

By  J.-F.  Ripoll,  A.  A.  Wray 


1.  Motivation  and  objectives 

The  two  stream  method,  also  called  the  two-flux  approximation,  distinguishes  between 
incoming  and  outgoing  radiation  for  improved  accuracy  relative  to  one-stream  methods 
while  retaining  reduced  numerical  cost  relative  to  full  RTE  solutions  (Mihalas  & Mihalas 
1984;  Siegel  & Howell  2001).  It  assumes  diflFerent  radiation  from  each  source,  each  uni- 
form over  its  half-space,  which  constitutes  the  weak  point  of  the  method. 

This  concept  has  recently  been  used  by  Dubroca  & Klar  (2002),  where  they  derive  a 
unidimensional  macroscopic  moment  radiation  model  for  each  stream;  they  call  this  a 
half  moment  model.  The  pressure  closure  is  obtained  by  using  maximization  of  entropy 
(Minerbo  1978),  which  allows  them  to  avoid  any  isotropic  assumption  about  the  radiative 
intensities. 

In  this  paper  a new  three  dimensional  half-moment  model  for  radiative  transfer  is  pre- 
sented for  a gray  medium.  It  describes  the  evolution  of  the  zeroth  and  first  directional 
half  moments  of  the  radiative  intensity.  The  closure  is  provided,  similarly  to  Dubroca 
& Klar  (2002),  by  the  maximum  entropy  concept.  This  work  generalizes  that  model  to 
three  dimensions. 

The  splitting  of  the  direction  of  propagation  fi  into  two  pieces,  and  f2“ , in  Dubroca 

& Klar  (2002)  was  done  by  cutting  the  fi-space  in  a static  sense,  meaning  that  the  same 
definition  of  + and  — was  used  at  all  points  in  the  domain.  This  direction  splitting  is 
clearly  the  best,  at  least  the  most  intuitive  one,  for  unidimensional  problems,  but  this  is 
not  necessarily  true  for  multi-dimensional  problems.  As  a matter  of  fact,  and  in  contrast, 
the  splitting  is  here  done  djmamically  according  to  the  direction  of  the  total  radiative  flux 
at  each  point.  At  any  point  of  the  domain,  our  model  considers  that  the  radiative  flux 
defines  the  main  direction  of  propagation,  the  positive  direction  and  a negative  one, 
in  the  opposite  direction,  This  dynamic  way  of  splitting  the  domain  of  directions 
appears  to  be  a natural  one  for  multidimensional  problems. 

This  particular  choice  for  the  splitting  also  allows  the  pressure  model  to  be  anal3dically 
computed,  which  is  not  the  case  for  a static  definition  of  the  half  moments.  However,  it 
does  have  a very  unfortunate  consequence:  since  the  half  spaces  fl+  and  fl~  are  depen- 
dent on  the  radiative  flux,  they  become  variable  in  time  and  space.  The  integration  of 
the  radiative  transfer  equation  (RTE)  over  these  subspaces  is  then  complicated.  Never- 
theless, if  the  radiation  is  assumed  isotropic  in  the  plane  perpendicular  to  the  direction 
of  propagation,  which  we  believe  is  a reasonable  assumption  for  a two-direction  model, 
the  integration  of  the  RTE  over  these  spaces  can  be  done.  Unfortunately,  the  integration 
introduces  unclosed  border  terms  involving  the  intensity  in  this  plane.  A closure  will  then 
be  provided  below  by  a model  for  this  intensity. 

In  order  to  derive  the  radiative  pressure  tensor  and  then  to  be  able  to  close  the  system. 
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a further  assumption  is  necessary:  the  positive  component  of  the  flux  is  assumed  to  be 
nearly  parallel  to  the  total  flux,  the  negative  one  becoming  then  necessarily  anti-parallel. 
This  is  obviously  true  in  one  dimension,  but  not  in  general,  and  constitutes  our  main 
assumption  in  the  pressure  derivation.  It  holds  nevertheless  exactly  in  two  important 
limits:  at  radiative  equilibrium  and  for  strongly  anisotropic  radiation. 

The  model  presented  here  (the  derivation  being  done  in  Ripoll  & Wray  (2003)),  called 
the  model,  is  a hyperbolic  system  consisting  of  a total  of  eight  equations  in  three 
dimensions,  four  equations  for  each  direction.  Each  half  model  has  the  classical  form  of 
a macroscopic  moment  model  in  which  the  pressure  tensor  is  constructed  from  the  well- 
known  Eddington  tensor  with  a particular  Eddington  factor.  Moreover,  different  source 
and  border  terms  occur.  The  latter  introduce  couplings  between  the  macroscopic  and 
microscopic  quantities  and  between  the  + and  — streams,  through  the  intensity  in  the 
plane  perpendicular  to  the  flux. 

One  of  the  major  advantages  of  this  model  is  that  independent  incoming  and  outgoing 
boundary  conditions  are  allowed,  which  is  not  possible  with  full  moment  models.  More- 
over, the  flux  stays  limited  by  the  speed  of  hght  and  the  underlying  intensity,  which  can 
always  be  deduced  from  the  macroscopic  quantities.  The  flux  is  described  by  a Planck 
function  at  radiative  equilibrium  or  by  one  (or  two)  Dirac  function  (s)  in  one  (or  two) 
direction(s)  of  propagation  in  the  anisotropic  Umit(s).  Furthermore,  the  maximum  en- 
tropy closure,  which  has  been  often  applied  to  radiative  transfer  (see  for  instance  (Ripoll 
2004))  will  be  shown  as  a very  useful  and  powerful  concept  allowing  the  derivation  of 
new  accurate,  well-defined,  and  robust  models. 

The  main  theoretical  application  of  the  half  moment  model,  treated  in  this  paper,  is  its 
reduction  to  a full  moment  model,  called  M^,  for  the  particular  but  important  case  of  a 
hot,  opaque  source  radiating  in  a cold  transparent  (or  semi-transpaxent)  medium  for  very 
specific  applications,  such  as  stellar  interiors  or  atmospheres,  or  combustion  problems. 
This  model  consists  of  four  equations  and  is  derived  from  the  half  moment  model  with 
fairly  simple  arguments.  The  model  is  tested  on  a simple  test  case  for  different  values 
of  the  opacity  and  will  be  shown  to  give  very  good  results,  better  than  those  obtained 
from  either  the  Pi  or  Mi  closures.  For  all  problems  presented  in  this  paper,  the  solutions 
obtained  by  the  new  models  are  compared  with  those  obtained  by  using  a ray-tracing 
solver  of  the  RTE.  It  will  be  shown  numerically,  partially  here  and  mostly  in  Ripoll  & 
Wray  (2003),  that  the  Mj*"  model  constitutes  an  improvement  of  the  existing  closures 
and  may  be  particularly  useful  for  treating  radiation  in  stellar  interiors  or  atmospheres. 

The  problem  of  a modeled  solar  atmosphere,  in  which  the  opacity  is  roughly  approx- 
imated by  that  at  a wavelength  of  500nm  (Vernazza  et  al.  1981),  will  be  solved  and 
discussed.  The  radiation  field  is  obtained  from  the  model  at  reasonable  accuracy. 

The  structure  of  the  paper  is  as  follows.  In  section  2,  the  model  is  presented.  In 
section  3,  for  the  particular  case  of  a hot,  opaque  source  radiating  into  a cold  medium, 
the  half  moment  model  is  reduced  to  the  model.  In  section  4,  we  first  solve  a simple 
and  academic  problem  to  validate  the  models,  followed  by  a simplified  solar  atmosphere. 
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2.  The  Half-Moment  Model 

2.1.  The  radiative  transfer  equation 
The  radiative  transfer  equation  in  an  non-scattering,  emitting,  and  absorbing  gray  medium 
is  given  by 

-dtl  + n • V7  = aB{u,  T)  - al,  (2.1) 

c 

■where  the  intensity  I — I{t,  r,  fi,  u)  is  a function  of  the  time  t,  the  position  r,  the  direction 
of  propagation  ft  and  the  frequency  u.  The  Planck  radiative  intensity  B describes  the 
isotropic  emission  of  the  medium  at  frequency  v and  temperature  T.  Here  c is  the  velocity 
of  light  and  a,  the  spectral  absorption  coefficient  or  opacity,  is  assumed  to  be  independent 
of  u.  By  integrating  the  RTE  over  frequency  and  introducing  the  quantity  J{t,r,fl)  = 
I{t, r,  fl,v)dv  = {I t,  we  obtain 

iSt J + n • VJ  - — T'*  -aJ,  (2.2) 

C 47T 

where  the  constant  a is  given  by  a = (87r®fc‘*)/(15h^c^). 

2.2.  Derivation  of  the  Half-Moment  Model 

Two  half  spaces,  defined  by  the  direction  of  propagation  of  the  radiation,  are  introduced, 
splitting  the  domain  fl  into  2 disjoint  pieces,  fi"*"  and  such  that  f2  = U fl~, 
n ft~  = 0,  and 

= {fi  / Fr/Fr  • n > 0},  (2.3) 

ft- = {ft  / Fr/Fr  • n < 0},  (2.4) 

where  designates  the  radiative  flux  and  Fr  its  normj:.  are  then  two  unit  half 
spheres,  obtained  by  cutting  the  sphere  ft  by  the  plane  perpendicular  to  the  flux  vector 
Fill  Pj.  = / Ffl  ■ S7  = 0}. 

The  three  first  moments,  Er,  Fr,  and  P^,  respectively  the  radiative  energy,  flux 
vector  and  pressure  tensor,  of  the  radiative  intensity  I according  to  the  frequency  and 


the  direction  are  defined  by 

Pfl=-(J)f2  FR  = {ftJ)^,  Fr  — -^{ft®ftJ)Q,  (2.5) 

In  a similar  way,  the  half  moments  Er^,  Fj^,  and  are  defined  by 

Er^  = , F^  = (fiJ)jj±,  P^  = - (f2  ® f2  J)jj±  , (2.6) 

By  construction,  the  following  properties  hold 

Er  = Er'^ Er- , Ffl  = F^  + F^,  Pfl=:P^4-Pfl.  (2.7) 

Three  unit  vectors  of  propagation  n,  n"*",  n~  are  defined  according  to  the  direction  of 
the  radiative  flux  and  half  fluxes; 

ii  = Fr/Fr,  n+  = Fj/P^,  n-  = F^/P^,  (2.8) 


It  is  assumed  herein  that  the  positive  component  of  the  flux  can  be  taken  to  be  approx- 
imately parallel  to  the  flux  itself,  n=*=  ~ ±n  (hypothesis  HI).  This  assumption  is  true 

t We  will  denote  the  integration  of  a function  / over  the  variables  X,  Y,  Z as  {f)x,Y,z 
J for  the  sake  of  brevity,  the  norm  of  a vector  v will  be  always  denoted  as  v. 
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when  the  radiation  is  isotropic  or  strongly  anisotropic  but  is  less  accurate  in  intermediate 
cases.  As  a consequence  of  the  previous  definitions,  the  negative  component  of  the  flux 
then  becomes  anti-parallel  to  the  flux  F/j.  The  previous  relationship  among  the  radiative 
fluxes  can  then  be  rewritten  as 

F«  = F+n+  + F^n-  ~ F+n  - F^n.  (2.9) 

Integrating  the  radiative  transfer  equation  (2.2)  over  each  half  direction  , leads  to 


dtJd^+  [ QVJdfi  = 


IcaaT^  - caE%. 


(2.10) 


Multiplying  the  radiative  transfer  equation  (2.2)  by  ft  and  integrating  over  the  half 
directions  fl^  leads  to 


- [ dtflJ  dfl  + c [ n ® n • V J dfl  = - ctF^.  (2.11) 

c 7n±  4 

The  two  subspaces  have  been  chosen  in  order  to  partition  the  domain  in  a physically 
natural  way  and  to  allow  the  radiative  pressure  to  be  closed  (see  section  2.3).  But  the 
partitioning  causes  the  bounds  of  integration  over  to  depend  on  the  flux  and  hence 
on  x,y,z,t,  and  so  disallows  commuting  derivatives  and  integrals  over  ft*.  Nevertheless, 
it  has  been  possible  to  perform  these  integrals.  The  following  relationships  hold  provided 
that  the  radiation  is  assumed  to  be  isotropic  in  the  plane  perpendicular  to  the  direction 
of  the  flux  (hypothesis  H2),  for  example  if  described  or  modeled  by  a radiative  intensity 
Ij_  such  as  that  in  the  next  section. 


/ dtJdflcidtE^, 

(2.12) 

f ft  • VJdft  V • Fp  - ttJxV  ■ n*. 

(2.13) 

Jn± 

[ dtflJdfl  ~ atFt  - -Jxdtn^, 

(2.14) 

ft  ® ft  • VJ  dfl  ~ V • P^, 

(2.15) 

where  Jx  = These  relationships  are  derived  in  Ripoll  & Wray  (2003)  using  as- 

sumption (H2)  which  leads  to  canceling  of  two  of  the  border  terms  in  (2.12)  and  in 
(2.15)  and  to  a simple  form  of  (2.13)  and  (2.14).  Thus,  the  model  is  given  in  three 
dimension  by 


dtE^  -f  V • F|  - ttJxV  • n*  - ^caaT*  - caE^, 

-dtFt  - - JxStn*  -1-  cV  • Pt  = -caaT^n^  - aFt. 
c c 4 


(2.16) 

(2.17) 


The  closure  of  the  model  will  be  provided  by  models  for  the  radiative  pressure  and  J± 
in  the  next  section. 
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Macroscopic  models  of  radiative  transfer 
2.3.  Closure  of  the  Half-Moment  Model 
2.3.1.  A model  for  the  pressure 

A special  radiative  intensity  I*  with  the  following  form  is  chosen 


2hu^ 


in  which 


T*(t, r,n) 


T 


and  where  A and  B are  defined  on  fl"'’  and  ft  as  follows 


A+n+ 

on  n • fl  > 0 

B+ 

on  n • fl  > 0 

and  B = 

1 

1 

on  n • fl  < 0 

B- 

on  n • fl  < 0 

(2.18) 


(2.19) 


(2.20) 


This  intensity,  determined  by  the  so-called  maximum  entropy  closure,  provides  a way 
to  close  the  radiative  pressure  Pfl  by  approximating  it  with  PJj,  computed  from  I*.  This 
constitutes  the  third  and  last  assumption  (H3)  used  here  to  derive  the  Mj  model. 

If  the  scalar  B and  the  vector  A are  defined  from  the  two  constraints  Eji  = {I*)vn 
Fii  = (ftl*)^^  f^,  then  the  intensity  I*,  defined  in  (2.18),  maximizes  the  radiative  entropy 
under  these  constraints  (Minerbo  1978).  The  subdivision  into  and  fl~  does  not  change 
this  property,  and  as  a direct  consequence  of  (2.7)  we  have  that  the  restricted  intensities 
/|jj±  maximize  the  radiative  entropy  under  the  constraints  Er^  = (/*),^_n±  = 

(fi/*)^  . Unfortunately,  it  is  not  possible  to  obtain  a closed  form  for  the  pressure  using 

these  constraints  directly;  instead  we  use  HI  to  approximate  the  constraints  as  follows 

^ and  (/‘(n^)),,^*  ^ (n/*(±n)),,„±  , (2.21) 

The  abbreviated  notation  7*(n^)  is  used  to  denote  I*  as  defined  using  A in  (2.20),  and 
7*(±n)  denotes  7*  using  instead 

A = ~ A"'‘n  on  n • H > 0 ; A = A“n“  ~ — A”n  on  n • fi  < 0,  (2.22) 

which  corresponds  to  the  use  of  HI  in  the  definition  (2.20).  The  computation  of  A,  from 
now  approximated  by  (2.22),  and  B is  shown  in  done  in  Ripoll  & Wray  (2003).  Actually, 
the  assumption  HI  was  not  needed  before  this  point  of  the  derivation. 

The  radiative  pressure  tensor  P^  is  approximated  by  P^*,  computed  in  Ripoll  &:  Wray 
(2003),  and  is  written 


P%=  f il®ftJdft~  [ H(g)HJ‘(n=^)dn  = P^*  ~ / ft  ® flJ*  (±n)  dfl 

~ PJe(f±,Ffl±)  = ==  D«(f±)Ffl± 

where  has  the  form  of  the  well-known  Eddington  tensor  evaluated  for  an  anisotropy 
= F|/(cE;?*),  given  by 


DR(f^)  = Id  + ^ 


f±2 


(2.23) 
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in  which  the  Eddington  factor  x is 


X{f^) 


8/f^ 

1+6/±  + -^1  + 12/±-12/±2 


if  > 0.5  and  xif^)  = 1/3  elsewhere. 


(2.24) 

There  is  then  no  coupling  between  outgoing  and  incoming  radiation  in  the  pressure  model. 
The  isotropic  and  anisotropic  limits  {f  0 and  / — > 1)  of  the  model,  which  define  its 
range  of  validity  are  discussed  in  Ripoll  & Wray  (2003).  The  corresponding  range  of 
is  [1/2,1];  this  range  will  hold  numerically  as  well  if  the  limitation  on  the  Eddington 
factor  X is  enforced  as  in  (2.24),  but  not  otherwise.  The  Eddington  factor  is  hence  in 
total  agreement  with  the  domain  of  definition  of  /^.  Without  this  limitation  on  x,  /* 
lower  than  1/2  implies  a radiative  intensity  defined,  problematically,  in  terms  of  an  A 
with  a negative  norm. 


2.3.2.  A model  for  I±  and  J± 

The  radiative  intensity  I*,  which  has  been  used  to  derive  mean  absorption  coefficient 
(Ripoll  et  al.  2001),  can  be  used  as  well  to  model  the  border  terms  coming  from  the 
integration,  since  I*  is  isotropic  in  the  plane  perpendicular  to  the  flux  and  we  have 


2hu^ 


fhvT  ' 


-1 


(2.25) 


since  • O = 0 for  fi  € Px  and  B'^  is  given  in  Ripoll  k.  Wray  (2003).  is  a Planck 
function  evaluated  in  when  radiation  are  isotropic  and  vanishes  when  radiation  are 
anisotropic  in  the  direction  of  the  flux.  The  main  problem  of  this  model  is  that  the 
intensity  is  reconstructed  from  the  plus  and  minus  macroscopic  quantities.  There  is  then 
two  intensities  7j[^,  which  should  have  the  same  value  in  P±,  but  nothing  can  guaranteed 
it.  The  definition  of  7j[  in  Px  in  then  not  obvious  and  we  propose  simply  an  average  of 
the  two  values 


T*  + 


T r*  -“X  + -^X 


T*  + 


Jx  ^ ^ 


+ j: 


where  is  given  by 


± _ aT^  _ ar|‘*  6(1  - A^f 


= 


47rP±^  47T  3 _ 3A±  + A^"^ 

3 ^4  + ^-12/±''  + 12/±  + 1)^ 


IGtt 


aT, 


R 


/±(6/±  + 1 + ^_12/±2  + l2/±  + 1)  ’ 


(2.26) 

(2.27) 

(2.28) 


is  a function  of  and  which  decrease  from  their  isotropic  values  = 

aT^^/2TT  (note  that  = aT'^/A'K  at  equilibrium)  to  their  anisotropic  values  j* 

0.  Moreover,  it  can  be  seen  that  = ErJAtt  could  constitute  a rough  model  for  Jx- 
One  could  think  the  average  (2.26)  is  arbitrary  and  any  average  of  the  form  ^ 
(aJj[^‘*'  + bJ'^~)/{a  + b)  could  have  been  chosen.  Actually,  a = 1 and  5 = 1 is  the  only 
combination  allowing  the  model  (2.26)  to  get  the  exact  isotropic  value  J±  = En/{^n)  = 
(P^  +P^)/(47t).  (The  demonstration  is  obvious  using  (2.28)  with  /*  — > 1/2  in  (2.26)). 
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Finally,  the  8 equations  in  three  dimension  of  the  Mj  model  are 

dtE^  + V • F|  - ttJIV  ■ n±  = - caE^,  (2.29) 

^ JlStn±  + cV  • (Dfl(f±)  Er^)  = ^caaT^n^  - aF%,  (2.30) 

and  have  been  derived  using  the  three  assumptions  HI,  H2,  and  H3.  The  positive  and 
negative  quantities  are  then  only  coupled  by  the  model  for  given  in  (2.26)-(2.28), 
which  constitutes  the  last  approximation,  and  when  the  full  moment  are  reconstructed. 


3.  A new  3D  moment  model  for  a hot  opaque  medium  emitting  in  a cold 
medium 

Prom  the  half  moment  model,  a new  moment  model  is  derived  for  a particular,  but 
nevertheless  important,  case:  a hot  opaque  source  emitting  in  a cold  medium.  It  can 
be  seen  in  (2.29)-(2.30)  that  the  non-linearity  of  the  radiative  pressure  does  not  allow 
reducing  a closed  full  moment  model  from  the  sum  of  the  two  half-moment  models  in 
the  general  case.  As  a matter  of  fact,  this  sum  leads  to 

dtEn  + V ■FR  = caaT‘^ -cctEr,  (3.1) 

iStFfl  + cV  • {■Dr{{+)Er+  + BRif-)ER-)  = -aFR,  (3.2) 

where  the  border  terms  have  canceled  when  added,  as  did  the  half  emission  terms  of  the 
flux  equations,  and  where  a pressure  remains  which  is  expressed  in  terms  of  unclosed 
quantities. 

Prom  this  system,  it  is  however  possible  to  reconstruct  a full  closed  moment  model  in 
the  particular  case  of  a hot,  opaque  opaque  source  emitting  in  a cold,  transparent  (or 
semi-transparent)  medium.  We  will  discuss  at  the  end  of  this  section  how  the  source  and 
the  exterior  medium  must  be  characterized  with  regard  to  their  temperature  and  opacity. 

We  now  discuss  the  two  domains  in  this  particular  problem. 

1)  First,  inside  the  hot,  opaque  source  where  radiation  can  be  considered  as  nearly 
isotropic  (/^  ~ 0.5),  the  linear  Pi  limit  holds  for  both  half-moment  models  since  in  this 
limit 

Pr  ^ ^i^/^(D+F+  + D«F«)  = \{E^  + El)  = \Er  (3.3) 

The  condition  ~ 0.5  implies  / ~ 0.  We  will  assume  for  this  model  that  radiation 
can  be  considered  isotropic  for  / < 0.5.  It  should  be  noticed  that  radiation  is  usually 
considered  nearly  isotropic  for  / < 0.3  (the  P\  validity  domain);  we  assume  that  it  is 
possible  to  extend  this  range  to  / < 0.5  without  too  much  impact  on  the  solution  because, 
for  many  applications,  the  anisotropic  factor  / is  predominantly  either  close  to  0 or  1. 
The  latter  case  will  be  handled  next. 

2)  Within  the  cold,  transparent  medium,  a cold  equilibrium  is  assumed  for  the  negative 
half  moments  while  the  positive  ones  are  strongly  anisotropic  due  to  the  radiating  source, 
implying  >>  0.5. 

The  exterior  medium  is  assumed  to  be  cold  enough  such  that  El  « and  ||Efl.  “ ||  << 
(or  indifferently  Fr~  « Fr^^  Vi  = 1..3).  Then  Er  ~ and  Fr  ~ F^  lead 
to  / ~ /+  and  hence  f » 0.5.  Finally,  using  these  approximations 

Pr  = Br+E+  + Br-EI  ~ Dfi(/+)E+  ~ D+(/)Efl. 


(3.4) 
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It  must  be  noticed  now  that  the  Eddington  factor  D^(/),  given  in  (2.23)-(2.24)  is 
limited  by  construction  to  = 1/3  when  /■•■  < 0.5,  and  since  /"*■  ~ /,  this  leads  to  the 
Pi  closure  when  radiation  is  isotropic.  This  tensor  then  naturally  makes  the  transition 
from  JDj  = 1/3  to  Z?^(/).  Thus  the  same  full  moment  model  is  found  to  be  valid  inside 
the  hot,  opaque  medium  and  the  cold  one  and  is  able  to  establish  the  transition.  The 
radiation  field  is  hence  fully  described  by  this  model,  called  here  , which  is  written 


dtEn  + V • Ffl  = caaT'^  - cctEr, 
iStFfl  + cV  • (D+(/)Ffl)  = -crFfl, 

with 

and 


(3.5) 

(3.6) 


(3.7) 


X+(/) 


8/^ 

1 + 6/  + V1  + 12/-12/2 


if  / > 0.5  and  x{f)  = elsewhere. 


(3.8) 


This  four  equation  model  gives  equations  for  the  full  moments  and  for  a range  of  the 
anisotropy  / € [0, 1],  It  should  be  noticed  here  that  this  reduction  to  a full-moment  sys- 
tem is  only  possible  because  the  positive  flux  always  radiates  in  the  full  flux  direction, 
whatever  this  direction  is.  With  a classical  spUtting  of  the  flux,  where  the  positive  com- 
ponent of  the  flux  is  statically  defined,  such  a reduction  is  not  possible  since  incoming 
and  outgoing  radiation  cannot  be  distinguished. 

This  model  describes  radiation  emitted  by  a source  which  is  isotropic  inside  (Pi  must 
be  valid  inside)  and  radiating  into  a exterior  domain  in  a strongly  anisotropic  way.  This 
signifies  that  this  exterior  medium  is  not  radiating  much  compared  to  the  source,  but 
might  for  example  absorb  part  of  emitting  radiation  of  the  source  provided  that  it  does 
not  re-radiate  towards  the  source.  Moreover,  if  the  exterior  medium  is  too  opaque  or  too 
hot,  the  anisotropy  factor  will  stay  lower  than  0.5  and  the  model  will  simply  reduce  to 

Pi. 

The  domain  of  validity  of  is  then  quite  large.  It  will  be  valid  in  many  applications 
where  the  main  interest  is  to  compute  the  heat  loss  of  a source  by  radiation  in  a non- 
remitting medium,  in  particular  for  luminous  flames  and  fires  burning  in  ambient  air  and 
for  stars  radiating  into  their  thin  atmospheres. 

We  will  see  in  the  next  section  that  this  model  gives  results  very  close  to  those  obtained 
by  a RTE  ray-tracing  solver. 


4.  Numerical  Results 

4.1.  A black  hot  sphere  radiating  into  an  infinite  semi-transparent  cold  gas 

The  first  case  is  a hot  black  sphere  of  radius  r = 1 which  emits  at  an  intensity  /^  = 1 into 
an  infinite  absorbing  cold  gas  (T  = 0)  of  two  different  and  small  opacities  a = 0.1,0.01. 
In  terms  of  macroscopic  quantities,  we  have  inside  the  sphere:  E^  = 2-k,  F^  = ±7t, 
Er,  = 47T  and  F^,  =0. 

For  this  case  the  positive  flux  will  remain  positive  and  the  negative  quantities,  which 
vanish,  do  not  need  to  be  computed.  Moreover,  the  exact  value  of  Jx,  which  here  is 
simply  Jx  = 0,  is  used.  The  outgoing  boundary  conditions  are  = 0,  drF^  = 0 at 
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Figure  1.  Radiative  quantities:  energy,  flux  and  anisotropy.  Left:  cr  = 0.1.  Right:  a = 0.01 


Figure  2.  Comparison  between  5 models  for  cr  = 0.1.  Left:  Radiative  Energy.  Right:  Radiative 

Flux. 

r = 1.6,  similarly  drEn  = 0,  = 0 for  the  full  moment  models,  in  order  to  simulate 

an  infinite  domain. 

We  compare  in  Fig.  1,  the  steady  states  of  the  radiative  energy,  flux,  and  anisotropic 
factor  given  by  M^'  to  those  obtained  by  the  RTE  solver  for  a = 0.1, 0.01.  Results  are 
in  very  good  agreement  for  each  of  these  quantities. 

We  compare  in  Fig.  3-2,  both  radiative  energy  and  flux  obtained  by  five  models  for 
these  opacities:  RTE,  Mi,  P\  and  the  new  M^.  Results  obtained  by  the  RTE  and 

the  M^^  solvers  were  presented  in  the  two  previous  figures.  The  full  moment  methods 
Ml,  Pi  and  are  solved  from  inside  the  hot  sphere  since  no  boundary  conditions  can 
be  prescribed  at  r = 1,  contrary  to  the  RTE  and  the  Mj^^  solvers.  is  in  its  range 
of  applicability  and  gives  results  almost  as  good  as  those  of  Mj^^.  It  might  be  presumed 
that  the  small  difference  between  them  comes  from  the  difference  of  the  boundary  con- 
ditions. In  all  cases.  Pi  gives  mostly  wrong  results  (for  cr  = 0.01,  the  anisotropic  factor 
/ = Fr/Er  found  by  Pi  at  r = 1.0  is  around  5,  meaning  a speed  of  light  five  times 
overestimated).  Surprisingly,  Mi  gives  a good  computation  of  the  flux  for  all  cases,  but 
a rough  computation  of  the  energy. 
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Figure  3.  Comparison  between  5 models  for  a — 0.01.  Left:  Radiative  Energy.  Right: 

Radiative  Flux. 


These  results  valid  the  models  and  let  presume  that  they  will  provide  an  accurate 
description  of  radiative  fields  of  the  solar  atmosphere,  which  presents  the  same  main 
characteristics  of  this  simple  geometry. 

4.2.  Computation  of  the  radiation  field  of  the  solar  atmosphere 

The  temperature  and  opacity  profiles  of  the  sun  are  taken  from  Stix  (1989).  These 
values  are  obtained  from  model  C of  Vernazza  et  al.  (1981)  and  are  representative  of 
a quiet  sun.  The  temperature  and  opacity  profiles  are  truncated  before  the  transition 
layer  to  the  corona,  which  is  here  considered  as  a zero  density  zone.  Confronted  with  the 
difficulty  of  determining  mean  absorption  coefficients  for  the  solar  atmosphere,  which  are 
needed  to  model  the  opacity,  we  obtained  the  opacity  from  the  optical  depth  r at  500  nm 
{a  = drsoo/dr). 

The  interior  of  the  sun  is  assumed  to  be  opaque  enough  to  be  represented  by  a black 
body  emitting  at  6910it’,  which  is  the  temperature  at  25  km  below  the  surface.  Outgoing 
boundary  conditions  are  identical  to  the  previous  case.  In  Fig  4,  the  in  spherical 
coordinates  with  spherical  symmetry  is  solved  and  compared  to  the  RTE  ray  tracing 
solution  for  the  solar  atmosphere. 

The  radiative  temperature  computed  by  is  coincident  with  the  one  computed  by 
the  RTE  solver  close  to  the  sun.  Far  away,  the  model  finds  a different  plateau  value 
(between  4960R'  and  4938R')  from  the  RTE  solver  (between  4700K  and  4650ii').  At 
r = 6.97  X 10®  m,  the  error  in  the  radiative  temperature  from  the  model,  which  finds 
4950jFf  instead  of  4682A',  is  5.72%.  At  r = 6.98  x 10®  m,  the  error  is  6%.  The  radiative 
flux  computed  by  both  models  presents,  close  to  the  sun,  a small  spike  followed  by  a 
smooth  bump  close  to  the  surface  which  is  due  to  the  opacity  profile.  The  error  of  the 
flux  computed  by  M1+  is  maximal  at  the  top  of  the  spike  where  it  reaches  17%,  while  at 
the  top  of  the  bump  the  error  is  around  10%.  Outside  this  highly  variable  zone,  the  profile 
of  the  radiative  flux  is  computed  with  a precision  of  5%  by  the  model.  Moreover,  it 
should  be  noticed  that,  in  this  particular  example,  radiation  is  near  isotropic  (/  ~ 0.5), 
a condition  in  which  the  error  introduced  by  is  maximal,  but  it  remains  as  we  can 
see  very  admissible. 
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Figure  4.  Radiation  field  of  the  sun  atmosphere.  Left:  Radiative  and  matter  Temperatures. 

Right:  Radiative  Flux. 


5.  Conclusion 

A new  three  dimensional  half  moment  model,  called  has  been  derived  using  a 

dynamic  definition  of  the  half  quantities  and  the  maximum  entropy  closure.  The  direction 
space  is  cut  into  two  pieces  according  to  the  direction  of  propagation  of  the  total  flux, 
at  each  time  and  position.  This  splitting  allows  the  radiative  pressure  to  be  closed,  using 
the  standard  maximum  entropy  closure,  with  approximate  constraints,  but  introduces 
non-trivial  border  terms. 

By  assuming  that  radiation  is  isotropic  in  the  plane  perpendicular  to  the  flux,  P±, 
which  seems  admissible  for  a two-flux  model,  the  radiative  transfer  equation  has  been 
integrated  on  these  moving  half-directional  spaces.  The  integration  is  difficult  but  leads 
to  small  and  simple  expressions;  the  border  terms  account  for  the  rotation  needed  to  set 
the  directional  space  in  the  direction  of  the  flux.  The  value  of  the  intensity  in  P±  must 
be  known,  and  a model  for  this  quantity  has  been  proposed. 

An  immediate  application  of  this  work  has  been  presented:  for  the  particular  case  of  a 
three  dimensional  hot,  opaque  source  emitting  into  a non  re-emitting  medium,  the  half 
moment  model  has  been  reduced  to  a new  full  moment  model,  called  here  M^. 

The  model  is  here  tested  for  the  case  of  a hot  one-dimensional  Gaussian  source  radiating 
into  media  of  various  opacities.  Very  good  agreement  has  been  shown  between  the  RTE 
and  the  model  independently  of  the  opacity.  The  relevance  and  the  correctness  of 
the  closure  has  been  shown,  in  one  dimension,  for  a case  theoretically  similar  to  the  solar 
atmosphere. 

The  radiation  field  of  a solar  atmosphere  with  an  approximate  opacity  has  been  com- 
puted with  the  Ml  model.  The  results  of  are  not  as  good  as  those  presented  in  the 
previous  theoretical  problem:  the  error  is  always  around  5%  for  all  unknowns,  except  in  a 
small  transitional  zone  where  the  flux  error  reaches  its  maximal  value  of  17%.  The  solar 
problem  close  to  the  surface,  being  near  isotropic,  is  a case  for  which  the  accuracy  of  M^ 
should  be  lowest,  but  is  nevertheless  found  to  be  admissible. 
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